## This code replicates the tables and figures for the Algeria survey portion of
## Matthew Cebul and Sharan Grewal, "Military Conscription and Nonviolent Resistance,"
## Comparative Political Systems, forthcoming.

## Load required packages
library(stargazer)
library(mediation)

## Load Data
data <- read.csv("surveydata.csv")

## Fix control prime as baseline
data$civprimes <- factor(data$civprimes, levels=c("civ_ct", "civ_conscript", "civ_corr", "civ_past", "civ_RU", "civ_UN"))




#################
## Regressions ##
#################

## Table 2 (and Appendix, Table 8)

one <- lm(prot_future~civprimes+age+female+edu+inc+amazigh+rural+unemp+student+partic+supp_regime+supp_military+supp_gendarm+supp_police+supp_parli+system+dembest+sharia+goals_1+mil_ad+boutef, data=data)

two <- lm(repress~civprimes+age+female+edu+inc+amazigh+rural+unemp+student+partic+supp_regime+supp_military+supp_gendarm+supp_police+supp_parli+system+dembest+sharia+goals_1+mil_ad+boutef, data=data)

three <- lm(soldiers_wont~civprimes+age+female+edu+inc+amazigh+rural+unemp+student+partic+supp_regime+supp_military+supp_gendarm+supp_police+supp_parli+system+dembest+sharia+goals_1+mil_ad+boutef, data=data)

four <- lm(officer_soldier~civprimes+age+female+edu+inc+amazigh+rural+unemp+student+partic+supp_regime+supp_military+supp_gendarm+supp_police+supp_parli+system+dembest+sharia+goals_1+mil_ad+boutef, data=data)

stargazer(one, two, three, four, single.row=T, digits=2, dep.var.labels=c("Will Protest","Military will repress","Soldiers won't shoot","Officers may, soldiers won't"), covariate.labels=c("Conscription", "Corruption", "Past Repression","Russia","United Nations", "Age", "Female","Education","Income","Amazigh","Urban","Unemployed","Student", "Have Protested","Supp Regime", "Supp Military","Supp Gendarmerie","Supp Police","Supp Parli","Change System","Supp Democracy","Supp Sharia","Goal: 5th Term","Military ad","Post-Bouteflika"))



###############
## Mediation ##
###############

## First, Baron & Kenny (1986): Appendix, Table 9

one <- lm(prot_future~civprimes+age+female+edu+inc+amazigh+rural+unemp+student+partic+supp_regime+supp_military+supp_gendarm+supp_police+supp_parli+system+dembest+sharia+goals_1+mil_ad+boutef, data=data)

two <- lm(prot_future~civprimes+age+female+edu+inc+amazigh+rural+unemp+student+partic+supp_regime+supp_military+supp_gendarm+supp_police+supp_parli+system+dembest+sharia+goals_1+mil_ad+boutef+repress, data=data)

three <- lm(prot_future~civprimes+age+female+edu+inc+amazigh+rural+unemp+student+partic+supp_regime+supp_military+supp_gendarm+supp_police+supp_parli+system+dembest+sharia+goals_1+mil_ad+boutef+soldiers_wont, data=data)

four <- lm(prot_future~civprimes+age+female+edu+inc+amazigh+rural+unemp+student+partic+supp_regime+supp_military+supp_gendarm+supp_police+supp_parli+system+dembest+sharia+goals_1+mil_ad+boutef+officer_soldier, data=data)

stargazer(one, two, three, four, single.row=T, digits=2, covariate.labels=c("Conscription", "Corruption", "Past Repression","Russia","United Nations","Age", "Female","Education","Income","Amazigh","Urban","Unemployed","Student", "Have Protested","Supp Regime", "Supp Military","Supp Gendarmerie","Supp Police","Supp Parli","Change System","Supp Democracy","Supp Sharia","Goal: 5th Term","Military ad","Post-Bouteflika", "Military will repress", "Soldiers won't shoot", "Officers may, soldiers won't"))



## Second, Imai et al (2010): Appendix, Table 10

## Mechanism 1 (Military will repress):
one <-  lm(repress~age+female+edu+inc+amazigh+rural+unemp+student+partic+supp_regime+supp_military+supp_gendarm+supp_police+supp_parli+system+dembest+sharia+goals_1+mil_ad+boutef+conscript_prime+corr_prime+past_prime+RU_prime+UN_prime, data=data[!is.na(data$prot_future),])

two <- lm(prot_future~age+female+edu+inc+amazigh+rural+unemp+student+partic+supp_regime+supp_military+supp_gendarm+supp_police+supp_parli+system+dembest+sharia+goals_1+mil_ad+boutef+conscript_prime+corr_prime+past_prime+RU_prime+UN_prime+repress, data=data)

med <- mediate(one, two, sims=1000, treat="conscript_prime", mediator="repress")
summary(med) 


## Mechanism 2 (Soldiers won't shoot):
one <-  lm(soldiers_wont~age+female+edu+inc+amazigh+rural+unemp+student+partic+supp_regime+supp_military+supp_gendarm+supp_police+supp_parli+system+dembest+sharia+goals_1+mil_ad+boutef+conscript_prime+corr_prime+past_prime+RU_prime+UN_prime, data=data[!is.na(data$prot_future),])

two <- lm(prot_future~age+female+edu+inc+amazigh+rural+unemp+student+partic+supp_regime+supp_military+supp_gendarm+supp_police+supp_parli+system+dembest+sharia+goals_1+mil_ad+boutef+conscript_prime+corr_prime+past_prime+RU_prime+UN_prime+soldiers_wont, data=data)

med2 <- mediate(one, two, sims=1000, treat="conscript_prime", mediator="soldiers_wont")
summary(med2) 


## Mechanism 3 (Officers may, soldiers won't):
one <-  lm(officer_soldier~age+female+edu+inc+amazigh+rural+unemp+student+partic+supp_regime+supp_military+supp_gendarm+supp_police+supp_parli+system+dembest+sharia+goals_1+mil_ad+boutef+conscript_prime+corr_prime+past_prime+RU_prime+UN_prime, data=data[!is.na(data$prot_future),])

two <- lm(prot_future~age+female+edu+inc+amazigh+rural+unemp+student+partic+supp_regime+supp_military+supp_gendarm+supp_police+supp_parli+system+dembest+sharia+goals_1+mil_ad+boutef+conscript_prime+corr_prime+past_prime+RU_prime+UN_prime+officer_soldier, data=data)

med3 <- mediate(one, two, sims=1000, treat="conscript_prime", mediator="officer_soldier")
summary(med3) # p=0.038



##################
## Supplemental ##
##################

## Driven by personally knowing a conscript

one <- lm(prot_future~conscript_prime_know1+conscript_prime_know0+corr_prime+past_prime+RU_prime+UN_prime+age+female+edu+inc+amazigh+rural+unemp+student+partic+supp_regime+supp_military+supp_gendarm+supp_police+supp_parli+system+dembest+sharia+goals_1+mil_ad+boutef, data=data)

two <- lm(repress~conscript_prime_know1+conscript_prime_know0+corr_prime+past_prime+RU_prime+UN_prime+age+female+edu+inc+amazigh+rural+unemp+student+partic+supp_regime+supp_military+supp_gendarm+supp_police+supp_parli+system+dembest+sharia+goals_1+mil_ad+boutef, data=data)

## Doesn't affect fighting back

three <- lm(fight2~conscript_prime+corr_prime+past_prime+RU_prime+UN_prime+age+female+edu+inc+amazigh+rural+unemp+student+partic+supp_regime+supp_military+supp_gendarm+supp_police+supp_parli+system+dembest+sharia+goals_1+mil_ad+boutef, data=data)

## Doesn't affect professional
four <- lm(prof~conscript_prime+corr_prime+past_prime+RU_prime+UN_prime+age+female+edu+inc+amazigh+rural+unemp+student+partic+supp_regime+supp_military+supp_gendarm+supp_police+supp_parli+system+dembest+sharia+goals_1+mil_ad+boutef, data=data)


stargazer(one, two, three, four, single.row=T, digits=2, dep.var.labels=c("Will Protest","Military will repress","Fight back","Military is professional"), covariate.labels=c("Conscription: Know Conscript","Conscription: Don't Know","Conscription", "Corruption", "Past Repression","Russia","United Nations","Age", "Female","Education","Income","Amazigh","Rural","Unemployed","Student", "Have Protested","Supp Regime", "Supp Military","Supp Gendarmerie","Supp Police","Supp Parli","Change System","Supp Democracy","Supp Sharia","Goal: 5th Term","Military ad","Post-Bouteflika"))



##################
## Demographics ##
##################

## Appendix, Table 7:
stargazer(data, digits=2, covariate.labels=c("Conscription", "Corruption", "Past Repression","Russia","United Nations","Conscription-Don't Know", "Conscription-Know", "Age", "Female","Education","Income","Amazigh","Rural","Unemployed","Student", "Have Protested","Will Protest","Military Will Repress","Soldiers Won't","Officers may, Soldiers won't","Supp Regime", "Supp Military","Supp Gendarmerie","Supp Police","Supp Parli","Change System","Supp Democracy","Supp Sharia","Goal: 5th Term","Military ad","Post-Bouteflika","Fight Back","Military is Professional"))



##################
## Verification ##
##################

## Appendix, Figure 6:
barplot(table(round(data$time,0)), xlab="Minutes", xlim=c(0,75), main="Length of Time to Complete Survey (N=1113)")

## percent match (code from here: https://github.com/andrewflowers/survey-fraud/blob/master/r_scripts/percentmatch.R )
percentmatchR <- function(matrix){
  matrix <- t(matrix)
  cols <- ncol(matrix)
  rows <- nrow(matrix)
  pmatch <- data.frame(entry=seq(1, cols), match=NA)
  for (c in 1:cols){
    comp <- matrix == matrix[,c]
    pmatch$match[c] <- max(colSums(comp[,-c], na.rm=TRUE)/rows, na.rm=TRUE)
  }
  return(pmatch)
}

test <- percentmatchR(data)

barplot(table(round(test$match*100,1)), xlab="Highest percent match for each respondent", ylab="Number of Respondents", main="Testing for Duplicates and Near-Duplicates", ylim=c(0,300))
mean(test$match>0.85)
